%% =======================================================================
% Figure 3_2
% ========================================================================
clear all; clc; close all;
rng(11);
NN    = 2e6;

currentFile = mfilename( 'fullpath' );
[pathstr,~,~] = fileparts( currentFile );
addpath( fullfile( pathstr, '_figures' ) );
addpath( fullfile( pathstr, '_aux' ) );

tt    = 1/15;
tx    = 0.01;
tx1   = 0.06; 
teps  = 1.60;
tepsb = teps; 
mu    = 0;
theta = sqrt(1/tt)*randn(NN,1);

noise_save       = [0.002:0.01:1.501];
delta_save_noise = zeros(1,length(noise_save));
delta_save_alt   = zeros(1,length(noise_save));

for ii=1:length(noise_save)
    teps = tepsb*noise_save(ii);
    
    v    = tx1/(tt+tx1);
    vs   = tx/(tt+tx);
    ty    = vs^2;

    ei1   = sqrt(1/tx)*randn(NN,1);
    xi1   = theta+ei1;

    wx    = tx1/(tt+2*tx1+ty*teps);
    wy    = vs*teps/(tt+2*tx1+ty*teps);

    ei2   = sqrt(1/tx)*randn(NN,1);
    eps   = sqrt(1/teps)*randn(NN,1);
    xi2   = theta+ei2;
    y     = v.*theta+eps;

    fi1   = v.*xi1+(1-v).*mu;
    fi2   = wx.*xi1+wx.*xi2+wy.*y+(1-2*wx-wy).*mu;

    error = theta-fi2;
    rev   = fi2-fi1; 

    beta  = rev\error; 
    delta = y\error;   

delta_save_noise(ii) = delta;

end

ratio = tx1/tx;
tx    = 3/4*tx;
tx1   = ratio*tx;
for ii=1:length(noise_save)
    teps = tepsb*noise_save(ii);
    
   v    = tx1/(tt+tx1);
    vs   = tx/(tt+tx);
    ty    = vs^2;

    ei1   = sqrt(1/tx)*randn(NN,1);
    xi1   = theta+ei1;

    wx    = tx1/(tt+2*tx1+ty*teps);
    wy    = vs*teps/(tt+2*tx1+ty*teps);

    ei2   = sqrt(1/tx)*randn(NN,1);
    eps   = sqrt(1/teps)*randn(NN,1);
    xi2   = theta+ei2;
    y     = v.*theta+eps;

    fi1   = v.*xi1+(1-v).*mu;
    fi2   = wx.*xi1+wx.*xi2+wy.*y+(1-2*wx-wy).*mu;

    error = theta-fi2;
    rev   = fi2-fi1; 

    beta  = rev\error; 
    delta = y\error;   

delta_save_alt(ii) = delta;

end

noise_ratio = noise_save;

make_figure_4_pretty(noise_ratio,[delta_save_noise;delta_save_alt])
saveas(gcf,'_figures/figure_3_lhs','epsc')

